Condensate fluctuations in the dilute Bose gas 
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The fluctuations of a number of particles in the Bose- Einstein condensate are studied in the grand 
canonical ensemble with an effective single-mode Hamiltonian, which is derived from an assumption 
that the mode corresponding to the Bose-Einstein condensate does not asymptotically correlate with 
other modes. The fluctuations are evaluated in the dilute limit with a proposed simple method, 
which is beyond the mean-field approximation. The accuracy of the latter is estimated; it is shown 
that the mean-field scheme does not work for the single-mode Hamiltonian, while for the Hartree 
Hamiltonian it allows us to estimate the condensate fluctuations up to a numerical factor. As a 
hypothesis, a formula is proposed that relates the fluctuations in the canonical ensemble with that 
of the grand canonical one. 
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I. INTRODUCTION 

The observation of Bose-Einstein condensation (BEC) 
in trapped alkali-metal gases [1] has stimulated theoreti- 
cal and experimental studies of the basic problems related 
to this phenomenon (see the reviews in [2-4]). In partic- 
ular, the statistical properties of the condensate are of 
especial interest for interacting bosons. In the case of 
the ideal Bose gas, fluctuations of a number of particles 
in the Bose-Einstein condensate (SNq) = (N§) — (No) 2 
have been studied thoroughly in a box [5] and a harmonic 
trap [6] as well. However, in the case of the interacting 
gas, this problem is rather subtle and requires a more re- 
fined approach. So far there has been no uniform treat- 
ment of this problem in the literature [7-19]. It is not 
even clear what the value of the power 7 is in dependence 
of the fluctuations on the total number of particles, 



(<5A> 2 > ex N\ 



(1) 



in various Gibbs ensembles (compare, e.g., the different 
results of Rcfs. [10] and [11, 12]). Apart from general 
theoretical interest, condensate fluctuations can be mea- 
sured, in principle, by means of a scattering of series of 
short laser pulses [20]. 

In this paper we deal with the problem of fluctuations 
in the grand canonical ensemble and consider the global 
gauge symmetry U(X) of the Bose system to be broken. 
In principle, the properties of the dilute Bose gas in the 
grand ensemble are described quite well by the Bogoli- 
ubov model [21], with the condensate terms being treated 
as asymptotic c-numbers and with three- and four-boson 
terms being neglected in the Hamiltonian. However, the 
statistical properties of the Bose-Einstein condensate are 
beyond the Bogoliubov theory due to the c-number re- 
placement. Consequently, in order to tackle this prob- 
lem, we need other approximations or assumptions that 
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keep the quantum nature of the condensate operators. 
Assuming that the condensate mode and other modes 
can be considered as quasi-independent, we obtain the 
usual thermodynamic fluctuations of the condensate with 
7 = 1 in Eq. (1) and show that this result is consistent 
quite well with the Bogoliubov theory. 

We consider a dilute Bose gas interacting with a short- 
range pairwise potential [i.e., the potential that goes to 
zero at r — > 00 as l/r m (m > 3) or faster]. The method 
used in our paper can be called the approximation of the 
quasi-independent mode. This assumption is nothing else 
but a generalization of Bogoliubov's relation (3) for the 
operators do and a (see Sec. II below). As a result, we 
consider the condensate mode and the other ones to be 
uncorrelated. In the framework of this approximation, 
one can easily derive the effective single-mode Hamilto- 
nian. Our consideration is primarily dedicated to the 
homogeneous case, except for Sec. IV B, where fluctua- 
tions in the Hartree model are discussed briefly in the 
nonhomogeneous case. 

The paper is organized as follows. In Sec. II we review 
some important issues in the field of Bose-Einstein con- 
densation, which might not be well known for a reader. 
In the next section we derive the effective single-mode 
Hamiltonian, which is employed for studying conden- 
sate fluctuations in the grand canonical ensemble. In 
Sees. IV A and IV B the fluctuations are evaluated in the 
framework of the mean-field approximation and beyond 
it. Possible corrections to the value of the fluctuations 
are discussed in Sec. IV C. In the last section the main 
results are summarized. 



II. GENERAL REMARKS 

Let us recall some important issues in the considered 
problem, following Bogoliubov's paper [22] . A number of 
bosons in the Bose-Einstein condensate can be defined as 
the macroscopic eigenvalue No of the one-body density 
matrix (ft (r)^(r')) - that is, / dV (r')^(r))M r ') = 



2 



N (j) (r), J d 3 r |0 o (r)| 2 = 1, where N is proportional to 
the total number of bosons in the thermodynamic limit, 
while the other eigenvalues are proportional to the unit. 
The eigenfunction <^n is a one-body wave function of the 
Bose-Einstein condensate, obeying the Gross-Pitaevskii 
equation in the first approximation. One can easily in- 
troduce the condensate operators do and dj by expanding 
the Bose field operators in the complete set of eigen- 
functions of the one-body matrix: r/;(r) = do</>o(r) + 
Sjjo%0j( r ); tnen > tne operator N = djdo describes 
the number of particles in the condensate: (No) = No- 
For a homogeneous system the index j is associated with 
momentum p, </> p (r) = exp(ip-r)/\/V (V is the volume), 
and the condensate operators coincide with the Bose 
zero-momentum operators. Bogoliubov noticed [21, 22] 
that in this case the operators ao/W and a^/W should 
be very close to c-numbers, since they commute in the 
thermodynamic limit V — > oo, n = N/V = const. This 
implies that 

where ip is an arbitrary phase and no = Nq/V denotes 
the density of the Bose-Einstein condensate. However, 
we always have the constraint (do) = (dj) = due to 
global gauge invariance, which is equivalent to conser- 
vation of the total number of bosons (i.e., \H,N\ = 0). 
Hence, the ground state is not stable with respect to 
an infinitesimally small perturbation that breaks the 
gauge symmetry. For a correct mathematical treat- 
ment of the homogeneous Bose system with broken sym- 
metry, Bogoliubov proposed [22] to include the terms 
— ^(dJe JV + aoe~' lip )y/V in the Hamiltonian, where the 
parameter v > and v — > 0. Now the absolute 
minimum of the energy corresponds to the state with 
(o>o)/W = (a\)* /\/V = ^n~o~e lip , because the gain in en- 
ergy per particle, by Eq. (2), is equal to 2vy / n /n in the 
limit V — > oo due to Bogoliubov's terms. We stress that 
the limit v — > should be performed after the thermody- 
namic one. These two subsequent limits yield the well- 
defined order parameter ($(r)) = (ao)/VV = ^/no~e lip 
with fixed phase <p. Thus, due to Bogoliubov's infinites- 
imal terms, the values of the anomalous averages (like 
(ip), (i>i>), and so on) change drastically, and such aver- 
ages can be called, following Bogoliubov, quasiaverages 
(hereafter by the term "average" we mean "quasiaver- 
age"). At the same time, the values of the normal av- 
erages (like ('ift'ip) and so on) do not change. In par- 
ticular, this is valid for the one-body (ip^ (r)ip(r')) and 
two-body (ri)^ (r2)ip(r 2 )ip(r' 1 )} matrices. Hence, in 
the thermodynamic limit the eigenfunctions of the ma- 
trices do not change when the symmetry is broken. On 
the other hand, the anomalous averages have a transpar- 
ent physical interpretation: (V>(r)) = v^o<A)( r ) is the 
eigenfunction of the one-body matrix associated with the 



maximum eigenvalue No, and (tp(r)tp(r')) is that of the 
two-body matrix with the eigenvalue N (N — 1) [23]. 
Hence, with the concept of broken symmetry we obtain 
a simple method of evaluating these eigenfunctions. 

Equation (2) can be derived from Bogoliubov's 
principle of correlation weakening [22], which takes 
a particular form (V'(r) • ■ ■ ^>t (ri) ••• V>(r 2 ) ••• ) ~ 
(V>(r)}(- • • V>t(ri) • • -V>(r 2 ) • • • } when |r - ri| — > oo, ••• 
and |r — r 2 | — > oo. Indeed, using the expression of 
the condensate operator in the homogeneous system, 
d, = / d 3 rip(r)/VV, we obtain, in the thermodynamic 
limit, 

(^-^('i)-fe)-) 

= 1 J d 3 r (^r)...^t( ri )...^( r2 )...J> 

= _L J d 3 r (^r))(...^( ri )...V;(r 2 )...} 

=M(...# (Pl) ...^ (ra) ...). (3 ) 

One can indeed replace the condensate operators by the 
c-numbers, and this procedure, according to Bogoliubov, 
still yields exact values [22] of all thermodynamic quan- 
tities in the limit V — > oo. Note that the above relations 
can easily be extended to the inhomogeneous case. For 
example, in order to consider the harmonic trap, it is suf- 
ficient to replace the volume V by l/u> 3 and define the 
thermodynamic limit as [24] N — > oo, Nut 3 = const (here 
uj = 3/Lu x Lu y uj z is the average harmonic frequency). 

We stress once more that the replacement of d and 
dj by the c-numbers gives always the main asymptotics 
of the averages of these operators (if the Bose-Einstein 
condensate exists), and this method is sufficient for ob- 
taining an exact value of any thermodynamic quantity 
(see Sec. 7 of Bogoliubov's paper [22]). However, the 
main asymptotic value is canceled when evaluating the 
condensate fluctuations (5 No), and, thus, we have to take 
into consideration here the quantum nature of the con- 
densate operators. Nevertheless, Eq. (2) results in the 
constraint 7 < 2 for the parameter 7 in Eq. (1). The 
case 7 > 1 implies that bosons in the Bose-Einstein con- 
densate and beyond it are strongly correlated, and the 
extent of the correlation is nonthcrmodynamically large, 
since the thermodynamic fluctuations of any extensive 
observable are proportional to the total number of parti- 
cles. The results 7 = 4/3 and 7=1 have been obtained 
in Refs. [10] and [11, 12], respectively, in the Gibbs canon- 
ical ensemble. The results of Refs. [11, 12] are consistent 
with ours for the grand canonical ensemble, while those of 
Ref. [10] are at variance with ours (see Sec. IV C below). 

If the symmetry is broken (that is, if the ground state 
is a superposition of states with different numbers of 
particles), then it is rather difficult to define what the 
Gibbs canonical ensemble is. Indeed, within a variational 
scheme, we cannot use the restriction N — const directly 
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and have to impose the additional condition (N) = N, 
which is equivalent to introducing the Lagrange term 
— pN in the Hamiltonian; this is nothing else but us- 
ing the Gibbs grand canonical ensemble. So, a reliable 
treatment of the canonical ensemble can be done in the 
framework of the number-conserving scheme (for recent 
developments of the scheme, see Refs. [25, 26]). For that 
reason, in this paper we restrict ourselves to the grand 
ensemble. 



III. APPROXIMATION OF THE 
QUASI-INDEPENDENT MODE 

Let us describe the model and methods of calculation. 
In order to evaluate the mean square of the number fluc- 
tuations in the Bose-Einstein condensate, it sufficient to 
know the condensate density matrix p c . According to 
general rules of quantum mechanics, it can be obtained 
from the total density matrix p by taking the partial trace 
in the occupation number representation 

Pc = ...Tr... P = (" " r V \P\ ■■■ ti p •")' ( 4 ) 

p=£G ■■■n p --- 

where we denote | • • • n p ■ ■ ■ ) = Ilp^o \ n v) an d the m dex p 
is associated with momentum in the homogeneous case. 
So p c depends on the condensate operators aj and a 
only and acts only on states of the condensate - say, the 
Fock states \N ) and their superpositions. By analogy 
with Eq. (4), one can also define the density matrix for 
all noncondensate states p out = Tr p, which depends on 

N 

the Bose operators a p and a p with p ^ 0. Certainly, it 
is difficult to evaluate the density matrix of the Bose- 
Einstein condensate directly from Eq. (4), and one has 
to use additional assumptions. 

As discussed in Sec. II, Eq. (3), a particular case of 
Bogoliubov's principle of correlation weakening, is valid 
in the thermodynamic limit. This suggests that the state 
of the condensate mode is, in effect, independent of the 
other modes. It means that one can use the approxima- 
tion for the density matrix p ~ /5 c /5 ou t- The parameters 
of the matrices p c and p out are assumed to be related in 
a self-consistent manner. The factorization of the total 
density matrix implies that the condensate and noncon- 
densate subsystems are considered to be uncorrclatcd. 
Hence, we have (AB) ~ (A) (B) for arbitrary A and and 
B, depending on the condensate and noncondensate op- 
erators, respectively. This is very consistent with the 
approximation [27] 

Ab~(A)b + A(b)-(A){b). (5) 

Let us apply the approximation (5) for a homogeneous 
Bose system with the pairwise potential V(r). In this 
case, the Hamiltonian reads, in terms of the Bose fields 



operators, 

h= J d 3 rv;t (r) ^^Z!_^ (r) 

+ J d 3 rd 3 r'V(\r~r'\)^(r)^(r')^(r')4>(r). (6) 

According to the above method, one should separate ex- 
plicitly the condensate operators in the Hamiltonian by 
means of substitution ip(r) = a a /^/V + t?(r) [here we 
denote $(r) = X) P ^o cx p(*P ' r )/W] and employ the 
decomposition (5) to every operator containing the prod- 
uct of the condensate operators and d or #t . As one can 
see, applying Eq. (5) to the Hamiltonian (6), we come to 
the approximation 

H~H C + H out + const, (7) 

and, hence, the density matrix of the grand canonical 
ensemble p ~ exp(— H/T) is really factorized (in this 
paper we include the term — pN in the Hamiltonian). By 
virtue of this factorization, the condensate fluctuations 
arise in the grand canonical ensemble due to the particle 
and temperature bath rather than the energy excitations. 
Here we derive the effective grand canonical Hamiltonian 
of the Bose-Einstein condensate: 

H c = ^V ^ + A(al 2 + al) + 2(B-p)ala 

+2C(al + a )Wy (8) 

where Vo is the zero component of the Fourier transform 
of the pairwise interaction V(r), and we introduce by 
definition the coefficients 

A = J d 3 rV(r)(d(r)S(0)), 

B = ^0(^(0)^(0))+ J d 3 rV(r){&(r)$(0)), 

C = J d 3 rV(r)(S1(r)S(r)S(0)). (9) 

In the above expressions, the gauge symmetry is certainly 
assumed to be broken, and for the sake of simplicity, we 
put ip = in Bogoliubov's terms (see Sec. II); in this 
case, all the coefficients are real. In general, they depend 
on the density and temperature and will be evaluated 
below. It is worthwhile to note that A, B, and C tend 
to constant values in the thermodynamic limit. As to 
the Hamiltonian H out , it is given by Eq. (6) but with c- 
numbers instead of the condensate operators [28], so we 
arrive at the standard asymptotically exact Hamiltonian 
with Bogoliubov's replacement. The coefficients (9) can 
be calculated by means of the Hamiltonian H out ■ 

When calculating the coefficients in the dilute limit, 
one can employ the results of our previous papers [23, 29- 
31], based on an analysis of the two-body density matrix. 
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One can show the following, 
(i) The functions 

ip(r) = 1 + ip(r), ifp(r) = \/2cos(p • r) + ip p (r) (p ^ 0) 

(10) 

are eigenfunctions of the two-body density matrix, which 
belong to the continuous spectrum and describe the scat- 
tering of pairs of bosons in a medium of other bosons. 
The quantum number p can be associated with the rel- 
ative momentum of that scattering [<p(r) corresponds to 
zero momentum]. Their scattering parts ip(r) and tp p , re- 
spectively, obey the boundary conditions ip(r), ippi r ) ~ * 
at r — ► oo. The Fourier transforms of the scattering 
parts can be expressed in terms of the Bose operators: 

m = ^ P (k) = ^ <^v±> , 



no 



2n 



where n p — (a p a p ) stands for the average occupation 
number of bosons. 

(11) The following limiting relations are valid: 

lmWr) = VMr), lim i<p(r) = <pM(r), (12) 

p— >0 n^O 

where (p^ (r) obeys the two-body Schrodinger equation, 
describing the scattering of two bosons in a vacuum with 
zero relative momentum. Now, using relations (11) and 

(12) , one can obtain the coefficients (9) in the dilute limit 
n — > after a small amount of algebra [32] 

A ~ n o [U {O) (0)-Vol B~2(n-n )V , 

C ~ 2V^(n-n )[C/ (0) (0)-y ]. (13) 

Here C/(°)(0) = / d 3 r <p(°\r)V(r) = Antfa/m stands for 
the scattering amplitude and a is the scattering length. 

Now all averages of the operators aj and So and their 
products can be calculated in the grand ensemble by 
means of the matrix 

p c = cM-H c /T)/Z c , (14) 

with the grand partition function for the condensate 

Z c = Trcxp(-H c /T). (15) 

Nq 

Note that the approximation of the quasi-independent 
mode was introduced earlier in Ref. [33]. The authors 
considered the quadratic terms only in operators a and 
fit in Eq. (8) (see the definition in Sec. IV B below) and 
obtained the condensate ground state, which turns out 
to be the coherent squeezed state in this case. However, 
this approach leads to additional approximations, which 
change drastically the final expression for the condensate 
fluctuations (see below). In order to avoid the additional 
approximations, we need to keep the four-boson Hartree 
term in the condensate Hamiltonian (8). 

The model (8) yields the asymptotically exact value of 
the condensate density provided that the coefficients (9) 



are known. Indeed, it follows from Eq. (3) that the con- 
densate operators can be replaced by their mean values 
(ao/W) = (al/VV) = ^/tm, in order to evaluate the 
main asymptotics of (H c ) in the thermodynamic limit. 
Then the condensate density can be found by minimiza- 
tion of (H c ) with respect to the c-number parameter 
no . The condition of this minimum can also be obtained 
from the initial Hamiltonian (6), since Bogoliubov's re- 
placement is asymptotically exact [22] due to Eq. (3). 
Hence, the condensate density no can be considered as 
a variational parameter whose value is obtained by min- 
imization of the thermodynamic potential of the grand 
canonical ensemble. With the well-known expression for 
an infinitesimal change of the potential, <50 = (SH), we 
have <9£l/9n = (8H/dn ) = d(H c )/dno = 0. Here the 
second equality is due to the condensate operators being 
involved explicitly only in H c . The last equation reads 



H = V Q n Q + A + B + C/s/rni. 



(16) 



We come to the equation obtained by Bogoliubov and 
treated by him as asymptotically exact [see Eq. (7.16) 
of Ref. [22]], since its derivation is based only on 
the c-number replacement of the condensate opera- 
tors. It can be rewritten as fi — J d 3 r'V(\r — 

r'\)(^{r')^{r')^{r))/^, where ^(r) = ^ + d(r). 
This equation relates the equilibrium value of the con- 
densate density no with that of the chemical potential. 



IV. CONDENSATE FLUCTUATIONS 
A. Method of parameter differentiation 

We adopt a standard method of evaluating the fluc- 
tuations from the Hamiltonian (8). One can derive a 
useful relation by differentiation of the expression for the 
average (No) — Tr (p c No) with respect to the chemical 
potential. If the Hamiltonian H c had commuted with 
N 0} we would have (5Nq)/No = T \driQ / 'd 'u)t / 'no from 
Eqs. (14) and (15), but the noncommutativity leads to 
corrections we need to estimate. One can apply the iden- 
tity cSexpf = Jo dz cxp[zf] 5fexp[(l - z)f] (here 6 de- 
notes variation) to the operator T = —H c /T and use the 
first two terms of the expansion 



e zl N a e- zl =No + [r,N ]- + [r,lr,N ]]- + --- . (17) 
This yields 



T— 

dp 

^Tr(e f iVo) 
dp 



T^=Tr(e*No), 



dp 



-l! dzT < 



No + ^(a Q 2 -ab 



C 



V(&1 - a ) 



e L N< 



■•)■ 



•5 



In order to estimate the second term on the r.h.s. 
of the last equation, we note that Tr (ctJ 2 e r iVo) = 

[Tr (al 2 e T N )]* = Tr (iV e r ao) owin g to tne reality of the 
coefficients in the Hamiltonian (8). Here we derive 

Tr ((aj 2 - a 2 )e f A>o) = Tr (e f [a 2 , A> ]) - 2(a 2 )Z 

and, in the same manner, 

Tr ((aj - a )e f A> ) = (a )Z. 

Thus, we arrive at the asymptotic expression 



(SN 2 ) = ^ /0«o\ 
AT n \dfiJ T 



2A 



C 1 



+ 



(18) 



We stress that the chemical potential (16) is involved 
in implicit form in the coefficients A, B, and C, but it 
is clear from the above consideration that one should 
hold the coefficients constant when taking the derivation 
d/i/dno in Eq. (18). As a result, we obtain, from Eq. (16), 



n (dfi/dn ) T = V n - C/(2y/no). 



(19) 



In the thermodynamic limit, Eq. (16) is the exact 
relation [22] for the chemical potential, and it can be 
applied [23] at any densities and for arbitrary strong 
pairwise potential [34]. By contrast, Eq. (18) cannot 
be employed in the strong-coupling regime, for which 
Vo — > +oo, since the divergence appearing in the co- 
efficients (13) is not canceled in this expression (18) 
[35]. However, one can use Eqs. (18) and (19) in the 
weak-coupling case, for Eq. (17) is in fact the expansion 
in terms of the coupling constant. Nevertheless, sim- 
ple physical arguments (see Sec. V in the first paper of 
Ref. [29]) can help us to extend our results to the strong- 
coupling regime. Since the properties of dilute quantum 
gases are ruled by the scattering length, the final expres- 
sion for condensate fluctuations should depend on the 
pairwise potential V(r) through mediation of the scat- 
tering length in the strong-coupling case. From this ex- 
pression one can derive the formula for the weak-coupling 
regime by means of the Born series for the scattering am- 
plitude (length): U^(0) = Airffa/m = U + U 1 + ---, 
where Uo = V and Ux = ~(2n)- 3 f d 3 fc V k 2 /(2T k ) < 0, 
Tfe = H 2 k 2 /(2m), and V k is the Fourier transform of the 
pairwise interaction. Otherwise, the relation obtained in 
the weak-coupling case should involve some first terms of 
the Born series, but coefficients before the term Uo = Vq 
are the same as before U^(0) in the strong-coupling 
case. Thus, in the "weak-coupling" formulas the sub- 
stitution V -> t/ (0) (0) (and U u U 2 , . . . -> 0) should be 
made to obtain the "strong-coupling" formulas. Perform- 
ing this substitution in the coefficients (13) yields A — > 0, 
C — > 0, and Eqs. (18) and (19) result in a simple final 
answer: 



m 



Nq Airtfno a ' 



(20) 



This equation is valid for sufficiently small depletion of 
the condensate - that is, when T < T c . Note that the 
above consideration allows us to avoid [29, 30] the diver- 
gence U\ — > — oo arising in the standard pseudopotential 
approximation V k = 4ttH a/m. 

The result (20) is a direct consequence of the single- 
mode approximation (8). Indeed, when deriving the 
Hamiltonian (8), we neglect the correlations between 
the condensate and noncondensate particles; besides, the 
condensate depiction is small for the dilute Bose gas. 
Consequently, we can put 



(SN 2 ) ~ (6N 2 ) = T 



dN 



mV T 
4irh 2 a ' 



(21) 



Here the familiar expression is employed for the fluc- 
tuations of the total number of particles in the grand 
ensemble and the formula for the chemical potential, 
= 4irH 2 an/m, in leading order at the density, which 
follows from Eqs. (13) and (16). Note that this expres- 
sion is nothing else but the relationship between the com- 
pressibility xt = (dn/ 'dP)r /n and the fluctuations of the 
total number of particles, where P stands for the pres- 
sure. Indeed, with the help of the thermodynamic rela- 
tion (dn/dP)T = n(dn/dfi)T and the definition of the 
compressibility, it can be written in the form 



(SN 2 ) ^ (SN 2 ) 



No 



N 



Tn X i 



For the ideal Bose gas a — 0, and we come to an infinitely 
large value of the fluctuations because of the nonphysical 
behaviour of the compressibility of the ideal Bose gas [5, 
6]. Note that the fluctuations of the total number of 
particles in the grand canonical ensemble remain finite 
and approach zero as T — > when a is finite, however 
weak the interaction may be; this is also valid for the 
one- and two-dimensional Bose gases [19]. 



B. Mean- field calculations 

It is interesting to compare the result (20) with the 
mean-field calculations for the condensate fluctuations. 
By separating the c-number part yfz in the condensate 
operator a = y/z + a (hence, (No) = N = (a T a) + z), 
the fluctuation of the condensate in the grand ensemble 
can be represented in the form 

(£A> 2 ) = (x + 1/2) 2 + y 2 - 1/4 + 2z(x + 1/2 + y), (22) 

where the notations x — (a) a) and y = (a^a^) = (da) 
are introduced, and the decoupling (a^ 2 a 2 ) = 2x 2 +y 2 is 
employed in accordance with Wick's theorem. 

Let us study the Hartree grand canonical Hamiltonian 



= 2~V a ° ~ Ma ° a °' 



(23) 



which is the model Hamiltonian (8) in the case of A = 
B = C = 0. To evaluate the parameters x and y in 



G 



Eq. (22), which correspond to the Hamiltonian (23), one 
can use the Gibbs-Bogoliubov inequality 



(24) 



Here (• • • )o means the averages in the grand Gibbs en- 
semble with the Hamiltonian Ho, and Q and f^o are the 
grand thermodynamic potentials corresponding to the ar- 
bitrary Hamiltonians H and Hq, respectively. Now, the 
basic idea is to choose H = and 

Ho = \ (A,a t2 + A*a 2 + 2B tfaj , (25) 

with arbitrary parameters Ao and Bo, and minimize the 
r.h.s. of Eq. (24) with respect to them. As a result, we 
find the stationary values of Aq and Bq. Note that Bo 
is always real, and we can put A = Aq if all coefficients 
are real in the Hamiltonian H. By using Bogoliubov's 
transformation, one can find the values x and y with the 
Hamiltonian (25) 



x = 



y 



Bo ., e 



1 



(a*a)o = ~z~ coth -— - „ . 
x 1 2s 2T 2 



(aa)o 



•— coth — , 
2e 2T 



(26) 



where e = \f B 2 — \ Ao\ 2 . We notice that at zero temper- 
ature this method is nothing else but the approximation 
of the coherent squeezed state (see, e.g., Rcf. [36]) for the 
condensate mode. 

It is more convenient to deal with the variables x and 
y rather than A and B , since (H)o is easily expressed 
via x and y with Wick's theorem and d(flo — (Ho)o) — 
— A dy — B dx. Hence, we come to the minimum condi- 
tions 

d(H) /dy = Ao, d{H)o/dx = B , (27) 
which should be solved in conjunction with 

d(H) Q /dz = 0. (28) 

For the Hartrce Hamiltonian we have (-ffh)o = ^F^ 3 ^ 
- Axz + 2zy + z 2 )- n(z + x), and Eqs. (27) and (28) 



yield 



y/ z coth ^EmE 



W-y/ z 

l + y/z 



coth 



(29) 



where the asymptotic formula V z/V ~ fi is utilized. 
The limits x/z — > and y/z — > at F ^ oo follow from 
Eq. (2), since z/V ~ n - If these limits were not fulfilled, 
the separation of the condensate operator into quantum 
and c-number parts would have made no sense. At zero 



temperature Eqs. (29) lead to x ~ — y ~ z 1 / 3 /2 4 / 3 and 
x + 1/2 + y ~ 2- 5 / 3 z- 1 /3 ; and Eq. (22) yields 



A 



25/3^1/3- 



(30) 



This result was obtained for the Hartree model by an- 
other method in Rcfs. [4, 8]. We note that the asymp- 
totics (30) is in agreement with Eq. (20) at T = 0, 
which reproduces only the exact limit {6Nq)/No = for 
Ao — > oo. At nonzer o tem perature we obtain in the same 
manners ~ -y ~ ^/T/fi^/z/2 andx+l/2 + y~ T/(2/z), 
and, hence, 



«) 

An 



3T 
2 7*' 



(31) 



The replacement Vo — » 47r H 2 a/m, discussed in Sec. IV A, 
results in Eq. (20) but with the factor 3/2. Certainly, the 
previous relation (20) is valid, as based on the more gen- 
eral considerations than the mean-field expression (31). 
Note that the method of the previous section, applied 
to the Hartree Hamiltonian (23), also leads to the same 
correct result (20). 

The mean-field scheme works even worse in the case of 
the Hamiltonian (8) , for the subtle balance of the terms is 
broken in Eqs. (22), (27), and (28). As a result, the term 
Vo = Uo vanishes in the limit V — > oo, and the main con- 
tribution comes in the condensate fluctuations from the 
term U\. Note that the approximations used in Ref. [33] 
lead to the same effect. Thus, the mean-field scheme, ap- 
plied to the model Hamiltonian (8), is not consistent for 
calculating the fluctuations, because it does not repro- 
duce the correct answer (20). Nevertheless, for a qualita- 
tive estimation of the condensate fluctuations, it is quite 
right to make use of the mean-field scheme, which repro- 
duces the correct answer up to a factor of 3/2. Note that 
the Hartree Hamiltonian can be easily written down for 
the nonhomogeneous system (see Sec. 2.8 in Ref. [4]); the 
coefficients depend on the condensate wave function (the 
eigenfunction of the one-body density matrix), which is 
the solution of the Gross-Pitaevskii equation. So it is 
possible to estimate qualitatively the condensate fluctu- 
ations for the trapped system by means of of Eqs. (30) 
and (31). 



C. Next-to-leading order corrections for the 
fluctuations 

In this paper we study the condensate fluctuations 
within the single-mode Hamiltonian (8). As is stressed 
in Sec. Ill, it means that we neglect the correlations be- 
tween the numbers of the particles in the condensate and 
beyond it; that is, we put [(A A out ) - (A ) (N out )]/V ~ 
in the thermodynamic limit, where by definition A out = 
N — No = J2 P ^o ^p^p 1S the operator of particles beyond 
the Bose-Einstein condensate. In Ref. [10] an opposite 
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idea was accepted that (SNq) = (<5A 2 ut ) due to the re- 
striction A = const in the canonical ensemble, but at the 
same time, Bogoliubov's scheme, which does not conserve 
the number of particles, was used there. This approach 
leads to nonthcrmodynamical fluctuations with 7 = 4/3 
in Eq. (1). As discussed in Sec. II, the replacement of the 
condensate operators by the c-numbers implies that Bo- 
goliubov's terms are involved in the Hamiltonian. Such 
a procedure leads to [H, A] ^ and can change the ac- 
cepted relation (SNq) = (<5A 2 ut ). On the other hand, 
within the conserving scheme (when [H, A] — 0) the 



Bogoliubov transformation b p = 



relates 



the creation ap and destruction a p quasiparticle oper- 
ators with not the particle operator a p but with b p = 
apSj/v^o [25, 26, 37]. Hence, the variance of the oper- 
ator J2 P ^o ^p^p i s no l° n g er equal to the variance of the 
number of noncondensate bosons (<5A 2 ut ). For this rea- 
son, the approach [10] is implicitly based on the assump- 
tion that [(A 2 A 2 ut ) - (A A ou t} 2 ]/A 2 ~ (A 2 ut ) - (A ou t> 2 
in leading order. A justification of this assumption is 
needed, since it may occur that the main contribution 
to the l.h.s. comes from the term (N^ ut (N§ - N§))/N§, 
which may be proportional to V 4 ^ 3 . The question re- 
mains open what approximation is more correct, the ap- 
proximation of the quasi-independent mode or the non- 
conserving approximation in conjunction with (SNq) = 
(6N% ut ). Note that our result (20) does not contradict 
to that of the papers [11, 12], in which the fluctuations 
were investigated in the canonical ensemble within the 
number-conserving simple scheme. We stress that the 
Bogoliubov model, based on the c-number replacement, 
is consistent with any value of 7 < 2, so all the ap- 
proaches [10] and [11, 12] and ours do not contradict 
to the Bogoliubov's c-number replacement. Thus, we are 
able to use the singe-mode Hamiltonian (8) until a de- 
cisive answer has been given what is the value of the 
correlator [(A^A\, ut ) — (A )(A ou t)]/"L within the number- 
conserving scheme [38]. 

Let us formulate the hypothesis for the interacting 
Bose gas that in the framework of the number- conserving 
scheme the relation (N N out ) /V ~ (N N out )c/V should 
be fulfilled; here, (• • • ) and (• • • ) c stand for the averages 
in the grand canonical and canonical ensembles, respec- 
tively This hypothesis is reasonable, since transitions of 
bosons from the condensate state to the noncondensate 
ones and back occur in the whole volume, and thus the 
boundary conditions seem to be of no importance here. 
In addition, we have the relations (NqN) c — AnA and 
(A A) - A A = T(dN /dfi) T ; the latter can be derived 
by differentiating (A ) with respect to the chemical po- 
tential (within the number-conserving scheme we do not 



face the difficulties discussed in Sec. IV A). Thus, from 
the accepted hypothesis we obtain 



An 



T 



dn 



+ 



An 



(32) 



v n \ u/j, j T iv 
Here the derivative (d[i/dno)T is not related to the 
formula (19), which concerns the single-mode Hamilto- 
nian (8), because the averages in Eq. (32) correspond 
to the full Hamiltonian with the pairwise potential; one 
should calculate it from the thermodynamic expression 
for the chemical potential. In particular, in the temper- 
ature region n J7 (0) (0) < T < T c [here T c is the tran- 
sition temperature of the ideal Bose gas and (0) = 
A-Kh 2 a/m] we have 



1 - 6V^ 



(noa 3 )V2T 



n !7( )(0) 



where (n a 3 ) 1/2 T/[nnJ7 (0) (0)] < 1 (see, e.g., Ref. [39]). 
On the other hand, the number-conserving approach of 
Ref. [11] yields 



An 



2v / 2\/^(n a 3 ) 1 



/2 



T V 
n o U(°H0)) 



in that temperature region [40]. Thus, with the proposed 
equation (32) we obtain the relation 



N rc t/(°)(0) 



l + (6 + 2V2)V^' 



(n a 3 ) 1 /2r 



n o U(°)(0) 



,(33) 



which contains the correction term in comparison with 
Eq. (20). 



V. SUMMARY 

Starting from the approximation of the quasi- 
independent mode, we derive the single-mode Hamil- 
tonian (8) and obtain the condensate fluctuations (20) 
for the grand canonical ensemble in the dilute limit. 
This relation is derived beyond the mean-field approx- 
imation. The mean-field scheme, applied to the Hamil- 
tonian (8), leads to incorrect results. For the Hartree 
Hamiltonian (23), the mean- field approximation results 
in Eq. (31) at nonzero temperature, which differs from 
the correct relation (20) by a factor of 3/2. With the 
help of the proposed hypothesis and the estimations of 
Ref. [11], the ncxt-to-leading term is obtained for the 
condensate fluctuations (33). 

This work was supported by RFBR grant No. 01-02- 
17650. 
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